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Abstract 

Recent neutron diffraction studies have shown that in Cr2(W,Te)06 systems, which consist 
of bilayers with strong antiferromagnetic inter-bilayer coupling between Cr moments, the intra¬ 
bilayer coupling between the Cr moments can be tuned from ferro (for W) to antiferro (for Te). Ah 
initio density functional calculations provide a microscopic understanding of the magnetic structure 
but cannot explain the magnitude of the ordered Cr^'*' moments. In order to understand the 
reduction of the ordered moment caused by quantum spin fluctuations (QSF) we have studied the 
magnon dispersion and the sublattice magnetization using a two parameter quantum Heisenberg 
spin Hamiltonian with tunable intra-(j) and antiferromagnetic inter- (J) bilayer couplings. In 
the j = 0 limit the system reduces to decoupled quantum dimers. The magnon dispersion and 
normalized sublattice magnetization (NSM) have been calculated using non-linear spin wave theory 
up to second-order corrections (I/S"^) in spin S. The effect of J is quite dramatic when j is 
ferromagnetic, a large peak in magnon dispersion appears at the corner of the Brillouin zone 
for large J (realistic case) in addition to changing the quadratic dispersion to linear. Inelastic 
neutron scattering experiments can conhrm this theoretical prediction. The effect of J on NSM 
is simple when j is ferromagnetic, NSM decreases monotonically with J. In contrast when j is 
antiferromagnetic the effect of J is to suppress QSF effect until J ^ j but further increase in 
J gives rise to increasing QSF reduction of the NSM. These theoretical results are qualitatively 
consistent with experiment. 

TAGS numbers: 71.15.Mb, 75.10.Jm, 75.25.-j, 75.30.Et, 75.40.Mg, 75.50.Ee, 73.43.Nq 
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I. INTRODUCTION 


Quantum spin fluctuations (QSF) play an important role in the low temperature proper¬ 
ties of quantum antiferromagnets (QAF), particularly in systems with low spin and low di- 
mensiond^^ For example, the ordered moment or sublattice magnetization (M^) in a nearest- 
neighbor (NN) id-dimensional spin-1/2 Heisenberg quantum antiferromagnet (HQAF) is zero 
for id = 1 (no long range order), 0.3067/1^ for D = 2^’® and 0.423/i5 for D = 3 (isotropic 
couplings)®, the latter two values obtained in a leading order approximation (to be discussed 
later in the paper). Thus QSF decrease with increasing S and increasing D. The interplay of 
QSF and covalency induced reduction of Mg in QAFs has also been a subject of great inter¬ 
est in the past, particularly in the parent compound of high Tc superconductors (La 2 Cu 04 
where the Cu^’*' ions have S = 1/2, form a square lattice and interact with NN isotropic 
Heisenberg AF interaction).^ 

Recently Zhu et al® have studied the magnetic structure of an interesting class of layered 
magnetic systems containing Cr®’*' (d® conhguration - spin 3/2) where the magnetic ions 
are arranged in bilayers. The inter-bilayer coupling is strongly antiferromagnetic due to 
the presence of Cr®''‘-Cr®''‘ dimers with strong intra-dimer AF couplings. The intra-bilayer 
couplings however can be either AF or ferromagnetic (F) depending on whether the system 
contains Te or W (also Mo) atom.® In fact the average intra-bilayer exchange can be tuned 
from one limit to the other in Cr2Wi_a;Tea,06 by changing x from 0 to 1. Neutron powder 
diffraction (NPD) measurements have determined the ground state spin structure and the 
values of the sublattice magnetization. The Te compound consists of antiferromagnetic 
bilayers which are coupled antiferromagnetically (AF-AF). In contrast the W and Mo analogs 
consist of ferromagnetic bilayers coupled antiferromagnetically (F-AF). In addition to the 
magnetic ordering the NPD measurements also give the values of Mg. The values of Mg 
are reduced from their atomic spin value 3.0//^ for Cr®"*" [Mg = g^sS) assuming g = 2 and 
quenched orbital angular momentum. This reduction can be due to covalency where the Cr 
d orbitals hybridize with O p orbitals and due to QSF.^ Electronic structure calculations 
reveal information about the reduction due to covalency whereas QSF-caused reduction can 
be calculated using a Heisenberg spin Hamiltonian, which is the major focus of this paper. 

Ah initio electronic structure calculations using density functional theory (GGA and 
GGA-|-U)^®^^^ correctly reproduced the magnetic ordering in these three compounds.^’® Since 
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the orbital angular momentum is quenched for the conhguration (three electrons in the 
local t 2 g orbitals; S = 3/2) the magnetic moment comes from the spin. The estimated 
exchange parameters (inter-bilayer or inter-dimer exchange J and intra-bilayer NN exchange 
j) were reasonable in view of the limitations of GGA or GGA-I-U approximations.® However 
the calculated values of the sublattice magnetization (~ 2 . 8 / 1 ^) was very close to the ionic 
value (S.O/is) indicating a small covalent reduction of the ordered moment.® In contrast, the 
experimental values are reduced to ~ 2.3/1^.® One possible reason for this reduction is QSF.^ 
Such dramatic reduction in ordered moment has been seen in many quasi-two dimensional 
QAFs, a classic example being La 2 Gu 04 which consists of antiferromagnetic 2D square 
lattice of S' = 1/2, where QSF reduces the sublattice magnetization (M^) by ~40%.^ In 
the present systems the reduction should be smaller (at least by a factor of three) due to 
F = 3/2. 

To visualize the magnetic ordering and magnetic exchange in the system we will consider 
the ground state spin ordering in Gr 2 (M=Te,W, Mo)06 (see Fig. 1).^® One has two bilayers 
(perpendicular to the 2 ;-axis) in the tetragonal unit cell (a, a, c) and four Gr spins/unit cell. 
The experimental unit cell parameters are a = 4.545A and c = 8.995A for Gr 2 Te 06 and 
a = 4.583A and c = 8.853A for Gr2W06.® The distance between the inter-bilayer (NN) 
Gr atoms i.e. Grl and Gr3 is 5 ~ 3.00A ra c/3, whereas the distance between intra-bilayer 
NN Gr atoms (Grl and Gr2 or Gr3 and Gr4) is ~ 3.80A. One bilayer contains Grl and 
Gr2 spins and the other contains Gr3 and Gr4 spins. The inter-bilayer coupling comes 
through Grl-Gr3 and Gr2-Gr4 dimers, it is antiferromagnetic and its strength is denoted 
by J. The NN intra-bilayer coupling is between Gr3-Gr4 and Grl-Gr2 and its strength is 
denoted by j. The magnitude of j is considerably smaller than that of J in these systems 
which can therefore regarded as weakly interacting quantum dimers, J being the strength 
of intra-dimer antiferromagnetic coupling. In Gr2Mo06 and in Gr2W06 the inter-dimer 
(also the intra-bilayer) coupling is ferromagnetic making the bilayers ferromagnetic. On the 
other hand, in Gr2Te06 the intra-bilayer coupling is antiferromagnetic. In thermodynamic 
measurements the high temperature properties (for example peak in heat capacity) are 
determined primarily by the dimers i.e. the energy scale is set by the intra-dimer coupling 
strength J whereas for the low temperature properties below the antiferromagnetic transition 
temperature the intra-dimer coupling j is responsible for the long range order. Also it plays 
an important role in magnon dispersion and quantum spin fluctuations. Experimental values 
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of the couplings are: J = —2.9 meV and j 
and i = 0.12 meV for Cr2W06.^^ 



—0.4 meV for Cr2Te06 and J = —3.8 meV 



FIG. 1. (Color online)(a) Schematic of the bilayer crystal structure of Cr2M06 (M=Te, W, Mo). 
Each Cr^'*' bilayer is separated by a (Te, Mo, W) layer.®’® (b) Positions of four chromium spins in 
the tetragonal unit cell of dimensions (a,a,c) are shown. The coordinates of the spins are: Crl: 
(0,0, c/2 — 5/2), Cr2: (a/2, a/2,5/2), Cr3: (0, 0, c/2 + 5/2), and Cr4: (a/2, a/2, c — 5/2). Within 
each bilayer Crl (Cr3) spin is coupled with Cr2 (Cr4) spin (shown by dashed lines with coupling 
strength j). On the other hand, for the inter-bilayer NN coupling Crl and Cr3 are coupled within 
the same tetragonal cell (shown by a thick solid line with coupling strength J) whereas Cr2 is 
coupled to Cr4 in the unit cell below and Cr4 is coupled with Cr2 in the unit cell above. 


The present systems are somewhat peculiar. The inter-bilayer coupling J is AF and 
strong. The intra-bilayer coupling j is small and can be either F or AF. If we put j = 0, 
then the system consists of non-interacting quantum spin S' = 3/2 dimers and the ground 
state consists of a product of dimer singlet states and there is no long range order (LRO). 
If on the other hand, J = 0, the system consists of non-interacting bilayers. When j is AF, 
the system is similar to the cuprates^ described above but the QSF reduces Mg by only 13% 
because S = 3/2. When j is F, there is no QSF reduction of Mg. An important question is 
how J and j interfere with each other. To study this interesting problem we have calculated 
the magnetic excitations and QSF for the interacting Quantum Spin dimer problem using a 
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Heisenberg quantum Hamiltonian. 

We investigate the difference in magnon dispersion and QSF between Cr2X06 (X=Mo,W) 
and Cr2Te06 at zero temperature using a two parameter (J, j) spin 3/2 quantum Heisenberg 
Hamiltonian. Our calculations use the well known Dyson^®’^® and Maleev^^’^® representation 
which facilitates the calculations significantly compared to the Holstein-Primakov represen¬ 
tation^®. The paper is organized as follows: we dehne the Hamiltonian in Section H and 
discuss the formalism using non-linear spin-wave theory^’^® along with results for the AF- 
AF bilayers (Section HA) and F-AF bilayers (Section HB). In Section HI we present the 
results for the magnon dispersion and sublattice magnetization and in Section IV the basic 
differences between the ground state properties when the inter-bilayer coupling is turned on 
between the two different types of long range ordered bilayers are discussed. We conclude 
our paper with a brief conclusion in Section V. 


II. FORMALISM 

It is well-known that quantum fluctuations play a signihcant role in the magnetic prop¬ 
erties^’^ and phase diagram of the system at zero temperature.Here we investigate 
the role of quantum fluctuations on the stability of the Neel state by calculating the magnon 
spectrum and see how these excitations reduce the sublattice magnetization from its Neel 
state value. We hrst express the fluctuations around the classical long range ordered an¬ 
tiferromagnetic ground state in terms of the boson operators using the Dyson-Maleev rep¬ 
resentation. The quadratic term in boson operators corresponds to the linear spin-wave 
theory (LSWT), whereas the higher-order terms represent interactions between magnons. 
We keep terms up to second order in 1/S. From the renormalized magnon Green’s functions 
and self-energies we calculate the magnon energy dispersion and the normalized sublattice 
magnetization (NSM) up to and including terms of order 1/S^. 

A. AF-AF Bilayer 

As discussed in Section H there are four chromium atoms in the tetragonal unit cell 
(shown in Fig. 1). Within each bilayer the spins of Crl and Cr3 are coupled antoferromag- 
netically with the spins of Cr2 and Cr4 respectively through the intra-bilayer NN coupling. 
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On the other hand, the NN Crl and Cr3 spins belonging to different bilayers are coupled 
antiferromagnetically within the same unit cell whereas Cr2 spin is coupled antiferromag- 
netically to Cr4 spin in the unit cell below and Cr4 spin is coupled with Cr2 in the unit cell 
above. 

The Hamiltonian of the system with bilayers with antiferromagnetic intra- and inter¬ 
bilayer couplings i and J (j, J > 0) takes the form 




"=jEEFd'sr+ssr'S, 

^=1 (hi) 


(4)A 

jn 


N. 


( 1 ) 


, 0(10 q(3)B , rq(2)B ^(40 cj(4)A c{2)B . 

' ^ 2 _^ 2 _^ “T 2 ■ ^m-1 “T 

n=l i 

where n represents the n-th unit cell index along the 2 ;-direction and i,j are nearest-neighbor 
sites within the same bilayer inside a unit cell. For example, represents the spin of 

Crl at site i in the n-th unit cell whereas represents spin of Cr2 in the (n -|- l)-th unit 

cell. A and B are indices for sublattices A and B. Within the bilayer spins of Crl, Cr4 are 
in A sublattice and Cr2, Cr3 are in B sub lattice. This spin Hamiltonian is mapped onto an 
equivalent Hamiltonian of interacting bosons by expressing the spin operators in terms of 
bosonic creation and annihilation operators a^, a for “up” sites on sublattice A (and b for 
“down” sites on sublattice B) using the Dyson-Maleev representation (dropping the indices 
i and n for the sake of illustration) 




(Li 


s+. = b] 


alttitti 

{2S) J 

(2S)J 


, S-., = s/2Sal Sh = S- aU, 


, = V^bj, = -s + b\b,. 


Bj 


.jUj. 


(2a) 

(2b) 


Substituting Eqs. (2) into Eq. (1) we expand the Hamiltonian perturbatively in powers of 
I/S' as 

H = H_i + Ho + Hi + --- , (3) 


where Hm is of order 1/5'”^“^. Note that H_i is just a number representing the classical 
energy. We do not discuss it further as it is irrelevant for the quantum fluctuations. Hence 
the 1/S expansion will be performed around the unperturbed Hamiltonian Hq which is the 
zeroth order Hamiltonian in this sense. Relative to Hq the terms Hi and H 2 are first and 
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second order terms respectively which are: 


H_, = -2jNS\A + t]), 
Ho = 
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-'in 


(4a) 


(4b) 


(4c) 


In Eq. (4a), the parameter rj = J/ j and N = N^NyN^ is the total number of unit cells. Next 
the real space Hamiltonian is transformed to momentum space using the Fourier transfor¬ 
mation for each £-th spin: 


W 

“in = 


^/N 




and then we use the following canonical transformations 


= exp 


ikz,c 

' 2 4 


5) Ok, b^fl = exp 




(5) 

( 6 ) 


where 6 is the inter-dimer separation (Fig. 1). Finally, we diagonalize the quadratic part 
Hq by transforming the operators Ok and b^ to magnon operators ttk and /9k using the 
Bogoliubov (BG) transformations 


“k = ^kojc + "ik/3-k, ^-k = m^al -F /k/3-k- 


(7) 










The coefficients and are defined as 


4 = 

where 


1 + Ck 
2ew 


1/2 


, rrik = -sgn(7k) 


1 - 
2ek 


1/2 


= -xJk, Xk = sgn(7k) 


1 - ek 


1 + ek 


1/2 


7 ik = cos(/c^a/2) cos{kya/2), 72 k = exp{-ik^c/2), 7 k = 


ek = (1-171.7'''". 

The resulting expressions in terms of magnon operators are: 

H_, = -2jNS\4 + r]), 

Ho = 2jS{A + r])^ [fk [aj,ak + /^Ik/^-k) + (ck - 1) 


471 k + V72k 
4: + 7] 


Hi = 
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(9b) 
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(10b) 
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r(2) 


1234 
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+ y2;34Q^Ia2/^l3/?l4 + y2Y34/^-l/^-2a3Cl4 + 14i2.34/^-4/^-3/^-2/^-l 
The coefficient Ak is given by 

Cp + |7pl ~ 1 
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-Y 

N ^ 


Ck- 


(lOd) 


( 11 ) 


In Eq. (lOd) momenta ki, k 2 , k 3 , k 4 are abbreviated as 1, 2, 3, and 4. The first term 
(Eq. (10a)) is the classical Neel energy, the second term (Eq. 10b) represents the excitation 
energy of the magnons within LSWT of 0{S^) and the leading order quantum zero-point 
energy, the term proportional to (ck — 1), and the third term (Eq. (lOd)) comprises 1/S 
contribution to the Hamiltonian. The first term in Eq. (lOd) is obtained by normal ordering 
the products of four boson operators with respect to creation and annihilation in the magnon 
operators, i.e. magnon creation operators appear always to the left of magnon annihilation 
operators. We have omitted the corrections to the ground state energy in Eq. (lOd) as they 
are not relevant to our calculations for the energy dispersion and magnetization. 

The quasiparticle energy «« Icir the a-magnon excitations measured in units of 2jS'(4-|- 
rf) up to second order in I/S' is given as 


= 2j^(4 + 7 ) 


Ck + 


:^k + 


(2.S) (2^) 


X 2 >(k. 6 k) 


( 12 ) 


9 























Expressions for the magnon Green’s functions and self-energies are given in Appendix A. 
The dynamic contributions to the second order self-energies are second order in the 
vertex factors These are the contributions which are missed by self-consistent spin- 

wave theory. Explicit forms of the vertex factors required for our calculations are given in 
Appendix B 1. 

Analytical expressions for the energy dispersion can be derived from the LSWT expression 
using Eq. (10b) for different regions in k space. They are: 


krtt 


EA/c^,0 ,0) = SjSsin { — ) \ r] + Acos'^ ^ T (0 < < - 






ky, 0) = AV2jS^/^iT2, {0<ky< -) 

a a 

-,k,)= aV2jS^/^, {0<k,< -) 


a a 


(13a) 

(13b) 

(13c) 


E^^{k,,ky,'^)=4V2jS^V + ‘^ 

E^{0,0,k^) = SjSy/rjsm(^^y {0 < k^ < ^) 


1 — cos^ 


,kra. 


cos 




TT, 


(0 < fca,, ky < -) (13d) 

CL 

(13e) 


The expressions for the sublattice magnetization are shown in Eq. (A7). In contrast to the 
quasi-particle energy, there is no hrst-order 1/S correction to the magnetization (Mi = 0 
for the AF-AF bilayer) as there are no off-diagonal i?k terms in Eq. (lOd). 


B. F-AF Bilayer 

In this case the NN intra-bilayer couplings j between Crl with Cr2 and Cr3 with Cr4 
are ferromagnetic whereas the inter-bilayer couplings J are antiferromagnetic between Crl 
with Cr3 within the same inter-bilayer and Cr2 spin with Cr4 spin between the adjacent up 
or down unit cells. 


N, 


^=1 {*j> 




Nz 


;(l)yl o(3)B 


rc(2)^ c(4)^ _1_ \ 

^in ‘ r)t^m ' ^in— 1 ‘ ^in ' ^m+lJ 


(4)S o(2)A 


(14) 


n=l i 

Here within the tetragonal unit cell Crl, Cr2 are in A sublattice and Cr3, Cr4 are in B 
sublattice. We then map the spin Hamiltonian onto an equivalent Hamiltonian of interacting 


10 










bosons using the DM transformation to obtain: 
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After Fourier transformation into the k space and using the canonical transformations 


(( 1 ) 
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(( 2 ) 

“k, “k = 6Xp 


ik^ c ■ 
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''‘1 = exp [-E(£ - 5 )] 6_^, jW = exp [E(£ _ i) 


“kj 

k) 


the Hamiltonian in Eq. (15) becomes: 

H_, = -2jNS\4 + r]), 

Ho = 2jS ^ I^KkCk (ak«k + /^Ik/^-k) + Kk(ek - 1) 

H, = 


J2§) ^ /^Ik/^-k) + 5k (alc/^lk + ttk/5-k) 
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1234 
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where 


7ik = cos{k^a/2) cos{kya/2), 72 k = cos{k^c/2), 

Kk = 4(1 - 7ik) + 7 , 7k = 

/^k 

Ck = (1 - |7k^^/^ 

. 2 1 ~ Ck ^ 7lk . 1 1 ~ Ck ~ 7k72k 

k k 

^k = - — [Ai{l - 7ik) + 7 ^ 2(1 - 7k72k)], 

Ck 

2 

Bk = —[Ai7k(l - 7ik) - 7^2(72k - 7k)]- 
Ck 


(18a) 

(18b) 

(18c) 

(18d) 

(18e) 

(18f) 


Explicit forms of the vertex factors E /234 are given in Appendix B 2. Once again, we 
have omitted the corrections to the ground state energy of 0(5'°) in Eq. (17d) as they are 
not relevant to our calculations for magnon dispersion and sublattice magnetization. 

The quasiparticle energy for a-magnon excitations, measured in units of 2jS up to 

second order in I/S' is given as 


= ^JS 


-I , 1 

Ek + TTTTTT^k + 


, Eli 


2Eu 


(19) 


( 2 ^) " ( 25)2 

Similar to the AF-AF bilayer, analytical expressions for the energy dispersion can be 
derived from the LSWT expression using Eq. (18c) for different regions. They are: 


E^(/c 3 ;, 0 , 0 ) = 8 jS'sin rj + AAr? (0 < < ^) 


E^^(-,/cj^,0) = 4\/2jS'^r/T2, (S><ky<-) 


pin 
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a a V21oV2/ c 
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, k„a , 


E {kx, ky, -) = 2jS 7 + 4|1 - cos(—) cos(-^) [> , (0 < k^, ky < -) 

k^c'' 


IT, 


E^(0,0,k^) = 2jSr)sm (0 < < ^) 


( 20 a) 

( 20 b) 

( 20 c) 

( 20 d) 

( 20 e) 


For the F-AF bilayer, the sublattice magnetization in Eq. (A7) has both first-order 1/5 
and second-order l/S '2 corrections. 
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III. RESULTS: ENERGY DISPERSION AND MAGNETIZATION 


A. Energy Dispersion 

We obtain the spin-wave energy dispersions for AF-AF and F-AF bilayers for S' = 3/2 
as a function of momenta {kx,ky,kz) for three different values of r; = J/j by evaluating 
Eqs. (12) and (19) in the hrst BZ. For the numerical summation we divide the hrst BZ in 
a mesh of points with Nl = 24. Then the contributions from points of p and q 
in the hrst BZ are summed up to evaluate the third terms in Eqs. (12) and (19). In the 
Dyson-Maleev formalism no cancellation of divergences occurs so that the convergence of 
the numerical results for ^ oo is very good.^"^ This is a crucial advantage over the use of 
the Holstein-Primakov representation. Our results for dispersion with Nl = 24: are accurate 
at least up to four signihcant digits. 

In Fig. 2 we show the energy-dispersions (for a-magnon excitations) calculated with 
LSWT (dashed lines), with 1/S corrections (dot-dashed lines in blue), and with l/S'^ cor¬ 
rections (solid lines) for the AF-AF and F-AF bilayers for rj = J/j = 5.0. In both cases 
(a) and (b) I/S' corrections are signihcantly large, but l/S'^ corrections are less than 0.5% 
compared to the results with 1/S corrections. It is interesting to point out that for the 
two-dimensional S' = 1/2 square lattice antiferromagnet with only NN AF interactions the 
dispersion at k = (vr, 0) and at k = (7r/2,7r/2) takes precisely the same value in linear 
spin-wave theory and with hrst-order 1/S corrections. However, in second order 1/S'^, a 
small dip appears which is about 1.4%. This number improves to 3.2% upon going to 
the third order 1/S'^.^^ Details can be found in Ref. (5) and (22). For our case the spin 
S being large, the higher-order corrections are signihcantly smaller and can be neglected. 
For the F-AF bilayer in Fig. 2(b) 1/S^ corrections are only noticeable between the region 
(vr/a, 0, 0) —(vr/a, tt/o, 0 ) where the dispersion is about 1% less than the dispersion calculated 
with only I/S' corrections. 

The spin-wave energy dispersions (for both LSWT and with 1/S'^ corrections) are shown 
in Fig. 3 for three different values of eta. For r; = 0 (J = 0) the bilayers are decoupled and 
we essentially have a decoupled 2D system. In this case we reproduce the existing results 
as shown in Fig. 3 - there is no energy dispersion between (0,0,7r/c) — (0,0,0). Note that 
for Fig. 3(a) the dispersion is shown in units of 2jS{4 + rj). So to compare the results with 


13 


1 




FIG. 2. (Color online) Spin-wave energy dispersion for the a-magnon excitations £^k,ao for the (a) 
AF-AF and (b) F-AF bilayers are plotted for rj = 5.0. Results from LSWT (dashed lines), with l/S 
correction (dot-dashed in blue), and with 1/5^ corrections (solid lines) are shown. For the AF-AF 
bilayer 1/S corrections harden the magnon dispersion and are significant (~ 10%) of LSWT results 
but 1/5^ corrections soften the spectrum but are hardly noticeable from the 1/S results (< 0.5% 
for most of the region except between (vr/a, 0,0) — {7r/a,7r/a, 0). For F-AF bilayer, 1/5 corrections 
harden the spectrum (by ~ 10%) from the LSWT results whereas the 1/5^ corrections lower the 
dispersion by about 1%. 

the analytical LSWT results for five regions given in Eqs. (13a) -(13e) we have to multiply 
our numerical data by the factor 2j5(4 -|- r]). On the other hand, the dispersion for the 
F-AF bilayer displayed in Fig. 3(b) is in units of 2jS. A peak at {7r/a,7r/a,7r/c) becomes 
more pronounced with increasing values of rj in this case. Similar to the AF -AF bilayer, the 
energy dispersions for the regions I-V can be checked against the analytical results obtained 
from LSWT as given in Eqs. (20a) -(20e). 

B. Sublattice Magnetization 

We calculate the normalized sublattice magnetization Ms/Mq from Eq. (A7) by numeri¬ 
cally evaluating Eqs. (A8a)-(A8c) as a function of rj. Especially to obtain the second order 
corrections M 2 we sum up the values of Nf/S points of k in one octant of the hrst BZ and 
Nf points of p and q in the hrst BZ, with Nl = 14 sites along one axis. The calculations are 
correct up to four signihcant digits. Fig. 4 shows the magnetizations (both LSWT and with 
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FIG. 3. (Color online) Spin-wave energy dispersion (for a-magnon excitations) E]^^aa for the (a) 
AF-AF and (b) F-AF bilayers are plotted for different values of 77 = J/j. Results from LSWT 
(dashed lines) and with 1/5^ corrections (solid lines) are both shown for ry = 0, 1.0, 5.0. 1/S 
corrections are significantly larger than LSWT results. However, 1/5^ corrections are not that 
significant compared to the results with 1/5 corrections. For rj = 0, the bilayers are decoupled 
and we reproduce the existing results for 2D Heisenberg antiferromagnet. For Fig. 3(b) a peak 
at {tt/ a^TT / a^TT /c) becomes more pronounced with increasing values of rj. For our plots we have 
chosen c/a = 2.0. 

1/5^ corrections) for AF-AF and F-AF bilayers as a function of ry = J/j. For Fig. 4(a) there 
are no first-order 1/5 corrections to the magnetization (as off-diagonal terms are not 
present in the quadratic part of Eq.(lOd)). Interlayer coupling initially enhances the magne¬ 
tization up to ry = 0.66 (from Ms/Mq = 1.3030 at ry = 0 to 1.3515 at ry = 0.66) and then with 
increase in spin-fluctuations magnetization decreases with increase in ry. 1/5^ terms in the 
Hamiltonian provide signihcant corrections to the LSWT results with increase in ry as the 
spin-fluctuations increase. Note that for ry = 0 spin-fluctuations are strong enough to reduce 
the scaled sublattice magnetization to 1.3030 from the classical value of 1.5 - similar to spin 
1/2 2D antiferromagnet on a square lattice where the sublattice magnetization reduces from 
0.5 to 0 . 3067.^’20 

For Fig. 4(b) as there is no antiferromagnetic inter-bilayer coupling for ry = 0 the bilayers 
are again decoupled and each bilayer is ferromagnetic. As a result there are no spin fluctua¬ 
tions. With increase in ry inter-bilayer AF coupling J introduces spin fluctuation corrections 
to magnetization which decreases monotonically. 1/5 terms in the Hamiltonian provide 
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significant corrections to the LSWT resnlts with increase in rj. Relative to I/S' corrections, 
corrections in this case are insignificant. 




FIG. 4. (Color online) Normalized sublattice magnetization M^/Mq for (a) AF-AF and (b) F-AF 
bilayers as a function of r/ = J/j obtained from LSWT (dashed lines), with 1/S (dot-dashed), and 
with corrections (solid lines). For Fig. 4(a) there are no first-order \/S corrections to the 

magnetization. Magnetization initially increases (spin-fluctuations decreases) up to 77 ~ 0.66 and 
then decreases monotonically as spin-fluctuations increase with rj. l/S^ terms in the Hamiltonian 
provide significant corrections to the LSWT results with increase in ry as the spin-fluctuations 
increase. In Fig. 4(b) when 77 = 0 there are no spin fluctuations as the 2D system is ferromagnetic 
(see text). With increase in 77 , the magnetization monotonically decreases. 1/S terms in the 
Hamiltonian provide significant corrections to the LSWT results with increase in 77 . 


IV. COMPARISON BETWEEN AF-AF AND F-AF BILAYERS 

We would like to discuss briefly the basic differences between the ground state properties 
of two different long range ordered bilayer states when the inter-bilayer coupling is turned 
on. When J = 0 (77 = 0) the bilayers are decoupled and the spin-wave dispersion is quadratic 
(linear) function of wave vector k for the F (AF) bilayers, and reduction of the normalized 
sublattice magnetization M^/Mq is 0 (~13%). When the antiferromagnetic inter-bilayer 
coupling J (77 > 0) is turned on, the changes in magnon dispersion and renormalization of 
the sublattice magnetization are dramatically different between the two cases. For the F 
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case magnons harden and the quadratic dispersion along the x|/-plane becomes linear for 
small k and for large rj the spectrum is essentially linear. A large peak develops around 
(vr/a, vr/a, tt/c) as the value of rj increases. In contrast for the AF bilayer where the magnon 
dispersion is already linear in the basal plane the introduction of AF inter-bilayer coupling 
does not do much. In fact it softens the magnon spectrum. The effect of non-zero rj on 
Ms/Mq is simple in the F case, Ms/Mq decreases monotonically with 77 , linearly for values of 
7] K, 2. For the AF case inter-bilayer coupling hrst increases Mg/Mo, as intra and inter-bilayer 
AF couplings compete with each other. When the inter-bilayer coupling becomes large the 
sublattice magnetization starts to decrease dramatically mimicking the weakly interacting 
quantum dimer limit. As an example, for the same value of ri{=5) the values of normalized 
magnetization Ms/Mq for the F(AF) case with corrections are 1.388 (1.236). 


V. CONCLUSIONS 

In this work we have investigated the interplay of inter-bilayer and intra-bilayer couplings 
on magnon dispersion and sublattice magnetization of system with S = 3/2. In particular, 
we studied higher-order quantum effects in a systematic perturbative spin-wave expansion 
in the inverse spin S. We have calculated the spin-wave energy and the magnetization up 
to and including the second-order corrections. They contribute signihcantly to the shape of 
the magnon spectrum and renormalization of the sublattice magnetization, especially as the 
ratio of interlayer and intralayer coupling increases. 
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Appendix A: Green’s fnnctions, Self-energies, and Magnetization 


The time-ordered magnon Green’s functions are defined as 

G„„(k,f) = -z(TK(t)4(0))), G;5/3(k,t) = -*(T(/3ljt)/3_k(0))), (Ala) 

G„/3(k,f) =-i(T(ak(t)/3-k(0))), G/3„(k, f) =-i(T(/3lj^(f)aj^(0))), (Alb) 


Considering Hq as the unperturbed Hamiltonian the Fourier transformed unperturbed prop¬ 
agators are 


G»^(k,a.) = Gj„(k,a.) = 0. 


1 

—oj — Ek + i5' 


(A2a) 

(A2b) 


with 5 0+. The spin-wave energy i?k = is measured in units of 2JS'(4 + rj) for the 

AF-AF bilayer and 2JS for the F-AF bilayer. The graphical representations of the Green 
functions are shown in Fig. 5(a). 

The full propagators Gjj(k, cj) satisfy the matrix Dyson equation 

Gjj(k, cg) = G°(k, cg) + ^G°^(k, cG)Smn(k, cG)G„j(k, cg), 

mn 

where the self-energy Sjj(k) can be expressed in powers of 1/(25) as 

S«(k,a.) = ^Elpk.!.,) + ^E<f (k.i.,) + .... 

The first-order self-energy terms read 

E«(k,a.) = E«(k.c) = ylk, (A5a) 

Ei«(k,a,) = E"(k,ia) = Bk. (A5b) 


(A3) 


(A4) 


The second-order self-energy terms originate from the Feynman diagrams in Figs. 5(b), 
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(a) ^ g P!! _ ^ _ _P 




FIG. 5. (a) The solid and the dashed lines correspond to the a and /3 propagators. Second-order 
diagrams for the self-energies Sao(k,w) and w) are shown in (b) and (c). The diagrams 

in (d) contribute only to SQ2(k, w). y(3)^ y(5)^ y(7)^ y(8) vertex factors (see main 

text). Note that at each vertex two arrows enter the vertex and two leave it which reflects the 
conservation of the total Sz component. 
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(c), and (d). The complete expressions read 
E<2(k,a.) = E®(k.c.) = (^) 


^( 3 ) 


1 \2 


N 


5^2(24Vk^[k+p-q])' 


Id 


( 2 ) 


Vr 


(3) 


k,p,q,[k+p-q] [k+p-q],q,p,k 


pq 


iu - Er, - En - E\ 


k+p-q] 


+ i5 


Vi 


(7) 


Vr 


(8) 


k,p,q, [k+p-q] [k+p-q] ,q,p,k 


CJ T Ep T/q T'jk+p—q] 7(5- 

E2](k,^^) = Eg](k..^) = (;i)'E2(2V,«|k+p-,l)" X 


(A6a) 


pq 


V 


( 2 ) 


V 


(7) 


k,p,q, [k+p-q] [k+p-q],q,p,k 


V 


(7) 


V 


(5) 


k,p,q, [k+p-q] [k+p-q] ,q,p,k 


(A6b) 

-OJ Ep T/q A'jk+p—q] T 7(5 UJ -[- Ep T/q T/[k+p—q] 7(5-1 
where [k + p — q] is meant to be mapped to (k + p — q) in the first BZ by an appropriate 
reciprocal vector G. In deriving Eqs. (A6a) and (A6b) we have used the symmetry properties 
of the vertices, see Eq. (B3a). 

The sublattice magnetization Mg for the A sublattice can be expressed as 


n/T — Q / t \ _ cx A C _L _L ^2 

Mg — S — {Ojj^QjjJ — D — l\o -|- “h 


(2^) (2A)2’ 


(A7) 


where 




_ 1 y- (k^k A 


Eu 


(A8a) 

(A8b) 




(*k+™k)^ 


^k777k ,^(2)/, p N 

“-^kj 


2(2/k/p/q/k+p-q) 


i‘i +™aC,.ik+p-,y 


( 8 ) 


k+p-q] [k+p-q],q,p,k 


pq 


(Ek + -Ep + Eq + Ek+p-q) 


Vr 


(5) 


^^k777kh^ p q [k+p-q] [k+p—q] ,q,p,k 

E^ - (Ep + Eq + Ek+p-q)2 


(A8c) 


The zeroth-order term AS* corresponds to the reduction of magnetization within LSWT, Mi 
term corresponds to the first-order 1/S correction, and M 2 is the second-order correction. 
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Appendix B: Vertex factors 


The vertex factors for the AF-AF and F-AF bilayers required for our calculations are 
given below. 


1. AF-AF Phase 


_ _ 
*^12;34 ~ ^ 


y(3) _ £ r 

'^12;34 — ^ 


yF) _ . 
*^12;34 ~ ^ 


i 

'^12;34 — ^ 


X37(2 - 3) - X47(2 - 4) - xiX 2 X 37 (l - 3) - 0 : 1 X 2 X 47(1 - 4) 
+ 2 : 1 X 27 ( 1 ) + 7 ( 2 ) + xiX 2 X 3 X 47(1 - 3 - 4) + X 3 X 47(2 - 3 - 4) , 

- Xi7(l - 3) - X27(2 - 3) - xiX3X47(1 - 4) - X2X3X47(2 - 4) 
+ xiX37(1) + X2X37(2) + xiX47(1 - 3 - 4) + X2X47(2 - 3 - 4) , 

- X2X3X47(1 - 3) - xiX3X47(2 - 3) - xi 7(2 - 4) - X27(l - 4) 
+ xiX47(2) + X2X47(1) + XiX37(2 - 3 - 4) + X2X37(1 - 3 - 4) , 

xiX47(1 - 3) + xiX37(1 - 4) X2X37(2 - 4) + X2X47(2 - 3) 

2 : 1 X 3 X 47 ( 1 ) - X2X3X47(2) - Xi 7 (l - 3 - 4) - X27(2 - 3 - 4) , 


'^12;34 — ^ 


xiX 47(2 - 4) + X2X47(1 - 4) + xiX37(2 - 3) + X2X37(1 - 3) 


- xi7(2) - X27(1) - xiX3X47(2 - 3 - 4) - X2X3X47(1 - 3 - 4) 


(Bla) 

(Bib) 

(Blc) 

(Bid) 

(Ble) 
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2. F-AF Phase 


y(2) _ 

'^12;34 — 


+ 


1 

2 L 


{ X 2 + a:iX 3 X 4 ){ 7 i (2 - 3 ) + 71 (2 - 4 ) + 71(1 - 3 ) + 71(1 - 4 )} 

{xi + a:2a:3X4){7i(l) + 71 (4 + 3 - 2 )} + (x2 + Xia:3a:4){7i(2) + 7 i (4 + 3 - 1 )} 
- 0:372(2 - 3 ) - 0:472(2 - 4 ) - 0:10:20:372(1 - 3 ) - 0:10:20:472(1 - 4 ) 


7 

4 L 


+ 


y(3) _ 

*^12;34 ~ 


+ 


0 : 10 : 272 ( 1 ) + 72 ( 2 ) + 0 : 10 : 20 : 30 : 472(4 + 3 - 1) + 0 : 30 : 472(4 + 3-2) 

^ - (0:3 + 0:1X20:4){71 (2 - 3) + 7i(2 - 4) + 7 i(l - 3) + 7 i(l - 4)} 
(x3 + xiX2X4){7i(1) + 7i(2) + 7i(4 + 3 - 1) + 71 (4 + 3-2)} 
0 : 33 : 472(2 - 3) + 72(2 - 4) + XiX2X3X472(1 - 3) + XiX272(1 - 4) 


y(5) _ 

*^12;34 ~ 


+ 


0 : 10 : 2 X 472 ( 1 ) - 0 : 472 ( 2 ) - xiX2 X372(4 + 3 - 1) - X372(4 + 3 - 2) 
1 


2 L 


- (X3 + xiX2X4){7i(2 - 3) + 7i(2 - 4) + 71(1 - 3) + 7 i(l - 4)} 
(x3 + xiX2X4){7i(1) + 7i(2) + 7i(4 + 3 - 1) + 71 (4 + 3-2)} 

0 : 10 : 30 : 472(2 - 3) - Xi72(2 - 4) - X2X3X472(1 - 3) - X272(1 - 4) 


+ 


^(7) _ 

'^12:34 ~ 


0 : 20 : 472 ( 1 ) + 0 : 1 X 472 ( 2 ) + X2X372(4 + 3 - 1) + xiX 372(4 + 3 - 2) 
^ (xiX 2 + X 3 X 4 ){ 7 i (2 - 3) + 7 i (2 - 4) + 71(1 - 3) + 71(1 - 4)} 
(xiX 2 + o:3X4){7i(1) + 7i(2) + 7 i (4 + 3 - 1) + 71 (4 + 3-2)} 
3;23:3{72(2 - 3) + 72(2 - 4)} + xiX3{72(1 - 3) + 72(1 - 4)} 


1 
4 L 


y(8) _ 

*^12;34 ~ 


0:1X3X472(1) - 0:20:30:472(2) - xi72(4 + 3 - 1 ) - X272(4 + 3 - 2 ) 

^ (xiX2 + X 3 X 4 ){ 7 i (2 - 3 ) + 7 i (2 - 4 ) + 71(1 - 3 ) + 71(1 - 4 )} 
(xiX2 + X 3 X 4 ){ 7 i( 1 ) + 7 i( 2 ) + 7 i (4 + 3 - 1 ) + 71 (4 + 3 - 2 )} 

I 0:1X372(2 - 3 ) + XiX472(2 - 4 ) + X2X372(1 - 3 ) + X 2 X 472(1 - 4 ) 
0:272(1) - 0:172(2) - X2X3X472(4 + 3 - 1 ) - XiX3 X472(4 + 3 - 2 ) 


(B2a) 


(B2b) 


(B2c) 


(B2d) 


(B2e) 
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Above, <I>G = exp(iG) where G is the reciprocal lattice vector. These vertex factors fulhll 
the following symmetry relations 


w(2) _ t/(2) 

y lO-'iA — ^ 19 - 


12;34 


12;43) 


w(3) _ T/(3) 


w(7) _ y(7) _ yy 

'^12;34 ~ '^21;34 ~ '^12;43i 


12;34 

h7) 


21;34i 


w(5) _ y 

'^12;34 — 


w(8) _ T/(8) _ Y 

'^12;34 ~ '^21;34 ~ '^12;43 


(5) 

21;34) 

( 8 ) 


If no reciprocal lattice vector is involved in the momentnm conservation, i.e., G 
are some additional symmetries 


(B3a) 
(B3b) 
0, there 


w(3) _ t^(5) . 

'^12;34 ~ '^12;34) 


y(7) _ t^(8) 

'^12;34 ~ '^12;34- 


(B3c) 
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